\documentclass{article}
\usepackage{Sweave}

\usepackage{graphicx,subfigure,amsmath,float}
\usepackage[margin=1in]{geometry}
\usepackage[dvipsnames]{xcolor}
\usepackage{pdfsync}
\usepackage{Sweave}

 
\SweaveOpts{echo=FALSE,fig=true,prefix.string=Figures/G,height=2,width=3}
\graphicspath{{Figures/}}
\setkeys{Gin}{width=0.45\textwidth}

\newcommand{\dottedline}{..}

\begin{document}

<<stuff, fig=false>>=
options(SweaveHooks=list(fig=function() par(mar=c(2.5,3,1,1), mgp=c(1.5, 0.5, 0), cex=0.7)))
source("../src/communitySim.R")
source("../src/simCancer.R")
source("../src/caseNumbers.R")
if(F) {
# how I got the 0.4899 figure
bob = function(qq) abs(qnorm(0.95, 0, sqrt(qq)) - log(10)/2)
optim(0.5, bob)
}

totalVar = 0.4899106
totalSD = sqrt(totalVar)
propGroupVar = 0.2


parameters = list(
	   geno=log(1.5), env=log(1.25), "geno:env"=log(2),
		 size=150000, probGeno=0.2, 
		 probEnv=0.3, genoMiss= 0.1,
		 envMiss = 0.1, oneEnvPerCommunity =F,
     cancerMissRate=0.15, falseCancerRate=0.00045, 
		 sdGroup= sqrt(propGroupVar) * totalSD, 
		 sdPerson=sqrt(1-propGroupVar)*totalSD, 
		 Enrollmentprobs = c(20000, 40000, 50000, 40000), 
		 agerange = c(35, 69),
		 Followup = c(5, 10, 20, 30),
		 Ncommunity =  50,
		 EqualCommunity=F
)
sig = 0.001

populationData = getPopData()  



deathFile=list(
          "F"=read.table("../data/MortalityFemale.txt",header=T),
          "M"=read.table("../data/MortalityMale.txt", header=T)
)


# common cancer means specific disease being selected. For cvd, we choose all of them  
CommonCancer=NULL

lambda =list(F=list(x=seq(30, 85, by=5), y= exp((-0.5)*0.585)*c(2.7, 3.9, 4.6, 6.6, 9.6, 12.3, 15.5, 18.2, 18.8, 
18.4, 16.4, 13.3)/1000),
M=list(x=seq(30,85, by=5), y=exp((-0.5)*0.585)*c(2.2, 3.7, 5.7, 8.5, 11.9, 15.8, 19.4, 23.3, 23.1, 22.3, 20.1,
14.7)/1000) )

DiffCancer=NULL
#DiffCancer= list(
#	"F"= matrix(cbind(lambda$F$y,lambda$F$y), nrow=length(lambda$F$x), ncol=2,
#dimnames=list(as.character(lambda$F$x),c("Total","diabetes"))),
#   "M"=matrix(cbind(lambda$M$y,lambda$M$y), nrow=length(lambda$M$x), ncol=2,
#dimnames=list(as.character(lambda$M$x),c("Total","diabetes")))
#	 )	
	 
propGroupVarseq = c(0.2, 0.35,  0.5, 0.65, 0.8,  0.95)
sdGroup= sqrt(propGroupVarseq) * totalSD 
sdPerson=sqrt(1-propGroupVarseq) * totalSD

SnumberCommunity = c(15, 30 ,50, 80)

Scolour = c("black", "black", "grey", "black")

thecolours =  apply(col2rgb(Scolour)/255, 2, paste, collapse=",")  
thelwd=c(1,1,2,2)
thelty = c(1,2,1,3)

thepch=1 

Nsim = 250

doSims=F
\section{Results}
\label{sec:results}

\subsection{Incidence Numbers}

<<caseNumbers,fig=false,results=tex>>=  
if(doSims) {
source("../src/caseNumbers.R")

          CommonCancer=NULL

caseNumbers = CaseCI(parameters, lambda, populationData, deathFile, DiffCancer, CommonCancer,verbose=F,
theGender=c("F","M"), Nsim, Percentile=c(0.025, 0.975)) 
  save(caseNumbers, file='diabNumCI.RData')

} else {
	load('diabNumCI.RData')
}
library(abind)
cases = abind(expected=caseNumbers$cases, caseNumbers$CI, along=1)
cases = aperm(cases, c(3,1,4,2))
cases =format(round(cases), width=4,justify='right')
cases=apply(cases, c(1,3,4), function(qq) paste(qq[1], '(', qq[2], ',', qq[3], ')',  collapse=''))
caseMat = matrix(cases, nrow=dim(cases)[[1]],
	dimnames = list(dimnames(cases)[[1]],
		NULL))
		
library(Hmisc)
caseMat = cbind(years=rownames(caseMat), caseMat)
latex(caseMat, file="",
	cgroup =c('Years', 'Colon', 'Stomach'), n.cgroup=c(1,2,2), 
	rowname=NULL, rowlabel=NULL,
#	rowlabel='years',
	colheads=c(' ', rep(dimnames(cases)[[2]],2)),
	caption='Average incidence numbers (with 95\\% prediction intervals) for males and females after various followup periods.',
	label='tab:caseNumbers', where='H', caption.loc='bottom'
)
@




\subsection{Graph}

<<genoEnvRR,fig=false>>=
if(doSims) {

 

genoEnvRR = seqPowerList(list("geno:env"=log(c(1.5, 1.75, 2, 3, 4))),
	parameters=parameters, CommonCancer=NULL, 
	SimulationTime=Nsim, lambda=lambda,
	populationData=populationData, deathFile) 
	
envRR = seqPowerList(list("env"=log(c(1.5, 1.75, 2, 3, 4))),
	parameters=parameters, CommonCancer=NULL, 
	SimulationTime=Nsim, lambda=lambda,
	populationData=populationData, deathFile) 
save(envRR, file="DiaenvRR.Rdata")	
cohortSize = seqPowerList(list("size"=1000*c(50,100,150,300)),
	parameters=parameters, CommonCancer=NULL, 
	SimulationTime=Nsim, lambda=lambda,
	populationData=populationData, deathFile) 
save(cohortSize, file="cohortSize.Rdata")	

save(cohortSize, genoEnvRR, envRR, file="DiaResults.RData")

} else {
load("Diacohortsize.RData")
load("DiaenvRR.Rdata")
load("DiagenoEnvRR.Rdata")
}
 

@
\begin{figure}[H]
\begin{center}


	\subfigure[Gene-Environment Interaction]{
<<allgenoenvRR>>=

toplot = genoEnvRR[,,'geno:env','All',as.character(sig)]



theLabels = rownames(toplot)
theLabels = gsub("[[:alnum:]|:|.]+=", "", theLabels)
theLabels = gsub(",", "", theLabels)
theLabels = signif(exp(as.numeric(theLabels)), 2)

 
  matplot(theLabels, toplot, col=Scolour, xlim=c(1,max(theLabels)),
  lwd=thelwd, lty=thelty,
  type="o",  xlab="relative risk", ylab="power", ylim=c(0,1),
  pch=thepch)

thecolours = toString(paste(dimnames(toplot)[[2]], 
c("( \\\\protect\\\\rule[3pt]{10pt}{1pt} )",
"( - - - )",
"( \\\\textcolor{Gray}{\\\\protect\\\\rule[3pt]{10pt}{2pt}} )",
"( \\\\ldots ) "
) )        )

@	
	}


\subfigure[Gene-Environment Interaction]{
<<Allgenoenvsize>>=
toplot = cohortSize[,,'geno:env',"All",as.character(sig)]  
 matplot(theLabels, toplot, col=Scolour, xlim=c(1,max(theLabels)),
  lwd=thelwd, lty=thelty,
  type="o",  xlab="Cohort Size", ylab="power", ylim=c(0,1),
  pch=thepch)  
@
}

	\subfigure[Environmental Effect]{
<<allenvsize>>=

@          toplot = cohortSize[,,'env',"All",as.character(sig)]
 matplot(theLabels, toplot, col=Scolour, xlim=c(min(theLabels),max(theLabels)),
  lwd=thelwd, lty=thelty,
  type="o",  xlab="Cohort Size", ylab="power", ylim=c(0,1),
  pch=thepc
}
	\subfigure[Environmental Effect]{
<<allenvRR>>=
toplot = envRR[,,'env',"All",as.character(sig)]
 matplot(theLabels, toplot, col=Scolour, xlim=c(min(theLabels),max(theLabels)),
  lwd=thelwd, lty=thelty,
  type="o",  xlab="Relative Risk", ylab="power", ylim=c(0,1),
  pch=thepch)
@
}

	\end{center}
\end{figure}



\end{document}
	